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Addendum for technology patent.txt 



Mr 393 



The. following public domain programs are used: 

1. BLASTN 2.0.3 [Nov-1 4-1 997 ] - from NCBI, NIH. 

2. perl, version 5.004_01, including the following includes from perl: 
Shell 

integer 
FileHandle 
strict 
Carp 

All the above are available through regular PERL distribution, and from 
the CPAN public .domain. 
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Clustering instructions 

Written by Raveh Gill-More, 26.2.98 

Updated by Avi Rosenberg, 28.5.98 



The clustering should be done with' the .user prod. 

The exec u tables are in /home/prod/bin / ARCH. ARCH is the 
architecture of the machine being' used - can be either DEC, 
SUN or SGI (but the last one isn't recommended, and many 
programs don't exist there). 
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(according to the current version 
je space (over 20GB recommended) . This 
>nnected to bioserv (in the rest of the 
3EDIR. In it make directories 
one, make directories named cluster and 



\ /...... . make a directory 



setenvGENBANKDIR/dir/bioserv 

c. Download to Israel (to the /dir/bioserv2/databases ) from Boston 

using the ISDN lines. Consult Artur about that. Shoud take -10 hours. 

2. UNCOMPRESS: Uncompress the- necessary files. You need the gbest* files, 
gbrod. seq.Z, gbpril.seq.Z and gbpr i 2 . seq . Z . 

cd $GENBANKDIR 

foreach f (gbest* gbpri* gbrod* ) 
echo $f 
gunzip $f 

end. 

3. FETCH: We need, to- fetch all the sequences of the organism we are interested 
in. This can be- done b.y running: 

fetch.pl organism=[ Ml H] molecule=(E|R] in= out = 

The first letter denotes the organism (Human or Mouse), the second is 
the s e qu e n c e - 1 y p e (Est or Rna), next comes the name of the 'GenBank file, 
followed by the output file name for the FastA format. 
You should work in two directories : "mouse" and "human". 
The commands are: . .1 
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A-3SS 
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t * 


. s 


eq 
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1 e 


c u 
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= E 



$BASEDIR/mouse/cluster/rna 
I sed ' s/ . seq//g ' ) 
:eq.out= $BASEDTR/mduse/clus 



end 



forhuman: 
cd$GENBANKDIR 

foreach i ( s \ls -1 gbpri* . seq I sed ' s/gbpri /'/g ' | sed 's/.seq//g' \ ) 

fetch.pl organism = H molecule=R in= gbpri${ i} .seq out= $BASEDIR/human/clus 

end 

foreach i ( 4 \ Is -1 gbest* .seq | sed *s/gbest//g'| sed *s'/.seq//g' x ) 

fetch.pl organism =H molecule=E in= gbest${i}.seq out=. $BASEDIR/human/clus 

end 

If you want to s pe ed up. 'the process, you can split up the fetching of the 
est files' to work' on several computers. 



MERGE: 
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e last 
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2 rna files ( one for mouse 


and. o 
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est file 
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using 


the script merge_and divide 


_f iles 



cd$BASEDIR/mouse/ cluster ' _ 

merge_and_divide_files.pl outpref ix=est max_size=25000 0000 file_type=fasta t 
cd$BASEDIR/human/c luster 

merge_and_divide_files.pl o u t p r e f i x = e s t m a x_ s i z e = 2 5 0 0 0 0 0 0 0 f i 1 e t y p e = f a s t a t 



Aft 


e r c h 


t m p 


.est* 


be 


delet 



5. DO_BLAST: Our sequences have to be scanned for repeats, vectors and highly 

expressed genes. This is. done using BLAST2. A script running BLAST and 
parsing the results is called " do_blast_cln.pl". Prepare the directory 

"blast_res" under the cluster directory of the human and mouse, and in it 
create the following directories: 

formouse: 

cd $ B AS E D I R / mou s e / c 1 u s t e r 
foreach' f (est? rna) 

foreach type (r v a) 

mkdir -p blast_res/${ f}_$type 

end 

end- 

for human: 

cd $ B A S E D I R / h urn a n / c 1 u s t e r 
foreach f (est? rna) 

foreachtype(rva) 

mkdir-pblast_res/${f)_$type 

end' 

end 



In each of these directories, one scan will be run. THe directory will end 
up containing one blast2 output file per query sequence (vector, repeat or 
abundantly occuring gene). Additionally, the files in' each directory will' 
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be parsed to produce the necessary annotation. 



Now, you need to format your data. Run: 

cd SBASEDIR /mouse/cluster 
.formatdb -i rna - p F 
f or each f .(est* ) 

f ormatdb - i $ f - p F 

end 

cd 3BASEDIR/ human/ cluster 
f o rmatdb-irna -p'F 
foreach 'f (est*) 

formatdb -i $f -p F 

end 
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REMARK: BLAST 2 is very sensitie to the abundancy of its query in the 
database. In some cases (such as ALU repeats) there is such a large number 
of hits that BLAST2 runs out of memory and crashes. This is the reason we 
keep the ests in a few separate files up to this stage. If we unite 
them before this stage - BLAST2 fails. Even so, BLAST2 takes a very long 
time and a lot of memory processing some of the more abundant repeats and 
genes. In the future it might be necessary to divide the human ESTs into 
more than two files. 

To check if any of the BLAST2 runs failed look at the length of' all" ".trap" 
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files which were created by the do_blast script. Sort them by length and 
see whether the shortest ones have terminated normally. All these files 
should be above IK bytes long if they are'OK.. The command is: 

ls-l*.tmp|sort-nr-k5|tail 

If any of the files is. incomplete - the BLAST would have to be re-run 
after dividing the appropriate EST file into two or more parts.. 

From now and on (until the assembly stage), you should work in the 
"cluster" directory of e a c horganism. 

6. MERGE: the pairs of EST bicfiles should now be merged while 
keeping the line's sort e d 'by EST name: 
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ileted. This is also the time 

ST files: 
. cat est? > est 

and remove est? (after checking the file sizes to make sure all is 
well). ' 

7. SORT_EST: The EST.s (mRNAs) should be sorted according to alphabetical 
order (to fit the order of lines in. the bicfile). To do that run: 

BigSort.pl in = rna out = rna. sorted -fasta keys=" ( ' A > (\ S+ ) \ s ' ) " 
BigSort.pl invest out = est. sorted -fasta keys=" ( ' A > (\ S + j \ s 1 ) " 

It is very advisable to run these commands on bioserv (i.e. locally on 
the machine that holds the files, and has a lot of memory) . 

Note : it is important to do the sorting to a new file and compare it's 
length to the original "est" file. Many things can go wrong while- sorting, 
and you don't want to re-create the "est" file. Only after you verify the 
' length of the two files are the same can you do: 

mv est. sorted', est' 

mv rna. sorted rna 

8. CLN: the cleaning stage. Assuming you have named all the files as in the 
previous sections you should do: 

cln di r= . in-rna -mrna 
clndir=.in=est 

All other parameters shouldn't be touched, but here is a list of what they 
mean: The suffixes control the names of the cleaning output files. 
The ," pbound" parameter is the bound for inserting annotation about repeats 
or- vectors to the sequence header, and the rest of the . numbers are 
parameters for the wraparound dynamic program which will try to find 
low-complexity regions in the sequence. It is very important to set the 
"-mrna" flag when cleaning mRNAs. Another parameter that' exists is the • 

"-no c u t _ n " flag, which - if specified - doesn't cut the ends of. the 

sequences based upon the number *of N's in "windows" of the sequence. 

Three output files are produced. They are (for the est file): 
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cat rna.xout.ok est.xout.ok > seq 

11. CLUSTER: Now.it is time to' cluster, or (actually) to contig. 

Use an alpha, with at least 256M of memory (currently mnm, godiva, tobler, 
bacio, panama and per u are good, but li.ndt and' chile are better, since they 
have more memory). You might want to change the " di r" and "in" parameters 
(if. they are not like the default' value - . d i r - . and i n = s e q ) , and the "mem" 
parameter. Mem should be set to 50 for human, 8 0. for mouse. Also you can 
play with the depth, and choose at most one of the flags (usually -sw, to 
perform a Smith-Waterman check at the end of the bin clustering). So, the 
command you will end up typing is: 



cluster mem = 50 -sw > & stdout . cluster 



The output is. a pairs file (seq. pairs) containing for each matching 
pair of sequences, data about the match (start and end in each sequence 
and quality) . 
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the 


oldest s 


e q u e 


n c 


e 


in the 


con 


t i 


g < 


a n 


d 




i 


f 


there are 


several s 


e q u e n 


ces 


from the 


same 


o 


1 d 


date 


- th 


e 


fir 


s t 




o f 




th 


e m 


lexi cogra 


p h i c a 


lly) 


, the da t 


e t h 


e 


s e 


quence 


was 


e 


n t e 


r e 


d 


i 


n 


t o 


GenBank, and 


the size 


of th 


e c o 


n t i g . Th i 


s f i 




w 


ill be 


upd 


a t 


ed 


1 a 


t 


e r 




t o 


contain more 


informati 


o n . 

































13. STATISTICS: we can now find out the sizes of our contigs: 

awk '{ print $2} 1 seq.clu | sort | uniq -c | awk '{print $1}' | sort +nr > s 



Now 


every 


1 


i 


n e 


in " s t a 


t " 


ho 


Ids 


t h 


e 


s i 


z e 


o 


f 


a si 


n 


gi 


e c 


o n t i g ( 


sort 


e d in 


d e s c 


ending 




o 


r d 


e r ) . See 


t 


h e 


stat 


i s 


t i 


c s 


a 


PP 


e n 


d i x 


a 


t 


the 


end f o 


r in 


structions 


on h 


o w t o 


o 


b 


t a 


in the r 


e 1 


e v a 


n t " 


c 1 


u s 


t e 


r i 


ng 




numb 


e 


r s 


f r 


om this 


f i 1 


e . 
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1 4 . 


GET_B I G 


_CONTI GS : 


The 


next 


few i 


terns are optional 


. I f 


t 


he biggest few 




c o n t i g s 


we found 


are 


VERY 


big , 


we would like to 


b r e a 


k 


them down further/ 




by r e - c 


lustering 


all 


the 


s e g u e n 


ces which partici 


pate 


i 


n th em, using 




more s t 


ringent c 


o n d i 


t i o n s 


. In t 


he .case of mouse 


( v e r 


1 


0 5) this was 




u n n e c e s 


s a r y but 


i n 


the c 


a s e of 


human it was. Th 


is d 


ep 


ends mainly on the 




c orap r eh 


ensivenes 


s of 


the 


a bu n d a 


nt gene database. 


So: 


i 


f you decide not 




to r e c o 


n t i g , you 


can 


go t 


o item 


19 ( R I C H _ F A S T A ) 


o t h e 


r w 


ise, we start by 




keeping 


a v e r s i o 


n of 


the 


first 


clustering output 









cp 


s e 


q • 


pa 


i r 


s 


seq. 


pa 


i 


cp 


s e 


q ■ 


cl 


u 


s 


eq.cl 


u . 


f 


cp 


s t 


a t 


s 


t a 


t 


. firs 


t 





Now we have to get the names of the seguences we want to extract. 
For human this i's: 

awk ' $ 4 > 500' seg.clu I awk '{print $1}' > seq. big .names 

For mouse use. 2 0 0 instead of 500. We now extract the sequences themselves 

get_ests.pl s eq .big. .names seq seq. big 



Rep 


e a 


t 


stage 


s 


y_f 


u 1 


1 _ 


o v e r 1 


a 


1 u" 


a 


nd 


" s t a 


t 



16. MERGE STATS: We need to merge several f.iles which are results of our two 
clustering runs. We begin by merging the " stat" files: 

awk '$1 <= 500' stat > tmp (200 for mouse) 

sort -m -nr trap stat.big >! stat 

Check everything is OK, and delete "stat.big" and "tmp". 

17. MERGE CLU: We now merge the "seq.clu" files: 

awk *$4 <= 500' seq.clu > tmp (200 for mouse) 

cat seq.big.clu tmp | sort -kb2 >! seq.clu 

The resulting file "seq.clu" will be sorted by contigs in alphabetical 
order. Delete "seq.big.clu" and "tmp", after checking all is well. 



1 8 . 


MERGE 


PA 




o v e r w r 


i t 




p r o d u c 


e . 




should 


b 



ra e 


r g e t h 


e 


" s e q . p a 


irs" f 


iles. 


Note: this will 


f 


the c 1 


u 


stering, 


which 


took 


several hours to 


ny 


doubt 




we might 


need 


the o 


riginal. output - it 


di 


f f e r e n 


t 


name be 


fore i 


t is 


overwritten ! 



delete_big_pairs.pl s e q . b i g . n a m e s seq. pairs seq. small .pairs 
cat seq. big. pairs seq. small, pairs >! seq. pa irs 

After checking, delete "seq. big. pairs" and "seq. small, pairs". At this 
stage you can also delete "seq" and "seq. big" . 

19. RICH FASTA: it is now (finally) time to create the rich fasta file 

that includes contig information though it will not be in its final form 
yet. Run: 

make_fasta.pl seq.clu rna.umsk.ok est.umsk.ok CCH1 0 6_ rich_fasta 
The parameter before last will act as a prefix to all contig names. 

20. CLONE INFORMATION: We now add information from the clone annotation of 
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our i 


n p u t 


sequences, t 


o u 


n 


i t e con 


tigs 


i n 


to clusters, w 


here we 


have 


no d i 


rect 


evidence of 


s e q 


u 


e n c e 


o v 


e r 1 a p 




but where at 1 


east t w i 


c e , lab 


t e c h n 


i c i a n 


s have assur 


ed 


u 


s t h a 


t 


a p a i 


r 


of sequences. 


one in e 


ach cluster, 


were 


p r o d u 


ced while se 


que 


n 


c i n g 


t h 


e s a m 


e 


cDNA. We run: 






c 


1 o n e_ 


clustersl .pi 


r i 


c 


h _ f a s 


t a 


1000 


2 


trap 






c 


1 o n e_ 


clusters2 . pi 


tm 


P 


c 1 o n 


e _ 


info 










t h e 1 


0 0 0 m 


eans we will 


no 


t 


a 1 1 e 


m p 


t t o 


u n 


ite clusters b 


igg er t h 


an 10 0 0- 


seque 


n c e s 


with other, c 


1 u s 


t 


e r s { 


t h 


e y t e 


nd 


to attract a 


lot of o 


ther clusters 


some 


of t h 


em erroneous 


i y ) 




The 


2 


means 


w 


e insist on 2 


inter — cl 


u s t <=* r 1 { n k ^ 


befor 


e we 


believe the 


inf 


o 


r m a t i 


o n 


( u n 1 


e s 


s o n e of th e c 


luster i 


s a 


s i n g 1 


eton, 


in which ca 


s e 


o 


ne li 


n k 


is a 


1 1 


wecan get an 


d it is 


considered 


to b e 


e n o u 


g h ) . t m p is 


a t 


em p o r a 


r y 


file 


name, which can 




t c u w 11 c (I 


the s 


e c o n d 


script has 


fin 


i 


shed. 


The re 


a s 


on the stage i 


s divide 


d i n t o 2 


scrip 


t s is 


technical ( 


mem 


o 


r y a 1 


1 o 


c a t i o 


n 


problems of pe 


r 1 ) . 




CLONE 


_ S T A T S : ' to obtain 


s t 


a 


t i s t i 


c s 


a b o u 


t 


what the clone 


i n f o r m a 


>- A U U ^ U 11 


tell 


us we 


run: 




















a 


w k ' { 


print $3}' c 


Ion 


e 


_i nf o 


1 


u n i q 




c | awk '{ prin 


t $ 1} • | 


sort + 0 n r > 


Now i 


n e a c 


h line. we ha 


v e 


a 


c 1 u s 


t e 


r siz 


e 


(number of con 


tigs con 


tigs). 


ADD CLONE 


INFORMATION : 


n o 


w 


we a 


dd 


c 1 us 


t e 


r information. 


from f i 


le 


clone 


_i n f o 


to seq.clu, 


a n 


d 


re-s 


o r 


t it 


by 


cluster name 


in alphabetical 



order : 



add_clone_info.pl seq.clu clone_info seq.clu. new 
sort - b -k 3,3 -k 2,2 -k 1,1 seq.clu. new >! seq.clu 



A f t e 


r a c h e 


c k , you can del 


e t e " s e 


q • 


clu.new" . 














The 


file no 


w has 7 fields 


per 1 i n 


e . 


They are: 


s e 


qu e n c e n 


a m 


e 


, c o n t 


i g 


name 


, c 1 u s t 


er name, sequen 


c e date 




contig size 


i 


n sequen 


c e 


s 


, c 1 u s 


t e r 


size 


in con 


tigs and cluste 


r size 


i n 


s e q u n e c e s . 














See 


the s t a 


tistics appendi 


x at t h 


e 


end for ins 


t r 


u c t i o n s 


o n 




how to 


obtain 


the 


r el e v a n 


t ".clone inform 


ation" 


n u 


rubers . 














SORT 


PAIRS: 


We now add the 


c o n t i g 


a 


nd cluster 


n a 


me to e v 


e r 


y 


line 


in 


"seq 


.pairs" 


, and sort itb 


y c 1 u s t 


e r 


name in al 


ph 


abetical 


o 


r 


d e r : 






s o r t _ p a 


irs.pl seq.clu 


seq . p a i 


r s 


seq. pairs. 


s o 


rted 










Then 


do * m v 


seq. pairs . sort 


e d seq. 


pa 


i r s * . 














COM P 


LE TE R I 


CH FA S TA : Now w 


e need 


t o 


add the cl 


u s 


t e r name 


t 


o 


the r 


i ch - 


fast 


a file. 


This is done 1 


ike thi 


s : 

















complete_rich_fasta.pl seq.clu rich_fasta rich^fasta . complete 



0 n 


t h 


e m o u 


s e 


da 


t 


a 


t h 


i 


s 




w 


o r k 


S 


fin 


e . 


H o w e v 


e r , 




o n 


t h 


e 




h 


u 


m 


a 


n 


d 


a 


t a 


the 


r e i s 


not 


e 


n o u g h 


m 


e m o 


r 


y 


t o 




r 


u 


n 


t h 


e 


s c r 


ipt 


in o 


n e 


g 


o . 


Th 


e 




t 


r 


i 


c 


k 


h 


e 


r e 


i s 


to divide 


the 


f 


i 1 e " 


s e 


q . c 


1 


u" 


i 


n 


t 


o 




two 


o 


r m 


ore 


part 


s . 


T 


o 


do 


t 


h 


a 


t 




t 


w 


o 


a 


dd i 


t i o 


n a 1 


par 


a m 


eters 


< 


fir 


s 


t 


a n 


d 




1 


a 


St ) 


m 


ay 


b e 


suppl 


ied 




t o 


th 


e 




s 


c 


r 


i 


P 


t , 




t el 


1 i n 


git 


w h i 


c h 


• line 


s 


o f 


»» 


s e 




c 


1 


u 




t 0 


u 


s e . 


So 


the 


who 


1 


e 


thi 


n 


g 




s 


h 


o 


u 


1 d 




be : 







complete_rich_fasta.pl seq.clu rich_fasta rich_fastal 0 500000 
complete_rich_fasta.pi seq.clu rich_f astal rich_f asta . complete 500001 

and then, after checking everything, do: 
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rmrich_fastal 

25. SORT: the Rich-FastA file needs to be sorted. We will use BigSort.pl 
again, but this time with, a different definition of the keys: 

BigSort.pl in = rich_fasta. complete out=rich_fasta -fasta keys=".( ' \ f CUVs+CC.-. 

It is very advisable to run this command on bioserv (i.e. locally on 
the machine that holds the files, and even more important, has a lot of 
memory). 

if the same size exactly', 



ich_f asta . co 


tnpl e te 


and rich_fasta 


a r 


are you can 


do: 






fasta . comple 


t e 






UT: we have 


two f i 1 


es which are p 


ass 


rich_fasta" 


(Rich-.F 


astA format, s 


o r t 


" seq. pairs" 


(also 


sorted, each 1 


i n e 


between two 


s e q u e n c 


es). The first 


s t 


the cluster 


to be 


in the same 'di 


r e c 


s done using 


the f 1 


ipper program: 





te assembly. They 

:ich_f asta" (Rich-FastA format, sorted by cl 

clones), and 

all of .the 

ly the right 

flipper gb_ver = your_version dir=. in = r i c h _ f a s t a p a i r s _ f i 1 e = seq. pairs ou 

27. ASSEMBLY: Since the assembly takes quite a lot of time, it is advisable 

to split up the. work between some roach i.n es. To do this, we pr e*p arealist 

of contigs that each machine will work on, taking into consideration the 

proportional strength of the machine (the DECs are about 2-2.5 times 
fasterthanSUNs): 

SplitClusters.pl in=seq.clu conf ig_f ile= assembly_config_file 

The file a s s e m b 1 y _ c o n f i g _ f i 1 e should be in the following format 

# any line beg ining with a # is a comment 
I 



# 


line 


s 


of the 


# 


type 


s 


of the 


#. 


line 


s 


of the 


t 


t h e i 


r 


types 


# 


for 


e x 


ample : 



DEC = 2 
S UN= 1 

# Note that the manchine types MUST be in upper-case 



b a 


c i o 


DEC 


go 


diva 


DEC 


m n 


m 


DEC 


1 i 


nd t 


DEC 


P e 


r u 


DEC 


pa 


n a m a 


DEC 


b i 


o s e r v 1 


SUN 


b i 


o s e r v 2 


SUN 


m i 


r 


SUN 



# m e a n s that DECs are twice as strong as SUNs , and we will use 

# 6 DEC machines (named bacio, godiva, mnm, lindt, peru, panama), 

# and. 3 SUN mach in es (named, bioservl, bioserv2, mir). 

The output of this are lists named .list. 
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Make a directory on each machine that you wish to run, preferably called 
. . . / p r odg b? ? ? / h u m a n / a s s emb 1 y / f i r s t _ r u n (for human), and 

. . .-/prodgb???/mouse/assembly/first_run (for mouse), if you don't have' 
a local disk on this machine, make a directory called 

$BASEDIR/hum'an/assembly//first_run (or use "mouse" for mouse), 
and cd to it. 

now open a window on each machine you want to run the assembly, and do: 
for human runs: 

setenv RUN -(usually run_name'is the machine name, but for 

instance if you have both bioservl and bioserv2 
that both run on bioserv, the name should 
reflect which run you are doing). 

setenv DIR$BASEDIR/human /cluster 

setenv' 'ORG human 

setenv GBVER ??? (accourding to the version) 
cp$DIR/$RUN. list. 

assembly dir= $DIR in= rich_fasta. flipped out=embl.$RUN transc_file=transcr 
for' mouse runs: 

setenv RUN . (usually run_name is, the machine name, but for 

instance if you have both bioservl and bioserv2 
that both run on bioserv, the name should 
reflect which run you are doing). 

setenv DIR $BASEDIR/ mouse/cluster 

setenv ORG mouse 

setenv GBVER ??? (accourding to the version) 
cp$DIR/$RUN.list. 

assembly dir= $DIR in= rich_fasta. flipped out=embl.$RUN transc_file=transcr 

In case any of the runs ends abnormally, you have to start it again from 
the. place it ended. To do this, you should look at the end of the "stdout" 
file of the run, and look for last appearance of something ' like this: 

############################################ #.################### 
### begin W 7 0 3 5 2 . A A 7 3 2 1 8 7 (size = 15) ### 

#**######**#*####*##########»########»########################## 

This means that the last contig the assembly tried to work on was called 
AA732187, in cluster W70352 (and it has 15 sequences in it). You should 
now make a new directory (for instance second_run), and rerun the assembly 
program, adding the parameter first-.. , v (where .... is the last contig 
that was worked on - in this case W 7 0 3 5 2 .. A A 7 3 2 1 8 7 ) . 

If this run also fails on this contig, look in the $RUN . list file, and 
see which contig comes after this, one, and try starting the assembly "run 
fromthisone. 



IMPORTANT NOTE: The assembly appends it's output to the existing files, 
so if you started an assembly run, and it was stopped (either by itself 
or manually), if you want to rerun it, you should either do it in a 
different directory (recommended), or delete the existing embl, cview, 
keywords, transcripts, graph and stdout files. 

It is very advisable to keep a list of the runs made in a file called 
Task. list at $ B A S E D I R / h u m a n / a s s e m b 1 y (or mouse), in the following format: 
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/A-V03 



# computer 

b a c i o 

godiva . 

m n m 

lindt 

peru 

panama . 

bioservl 

b i o s e r v 2 

m i r 

bac.io 



dir 

/dir/bacio2/prodgbl07/human/assembly/bacio/ first run 

/dir/godiva2/prodgbl07/huraan/assembly/godiva/fi rst run 

/dir/bioserv2/pr\odgbl07/human/assembly/mnm/first_run 

/ d i r / 1 i nd t 3 / p r odgb 1 0 7 / human / a s semb 1 y / 1 i nd t / f-i r-s t _r u n 

/ dir/perul/prodgbl07/human/assembly/peru/ first_run 

/dir/panamal/prodgbl07/human/assembly/panama/first run 

/dir/bioserv2/prodgbl07/human/assembly/bioservl/first_run 

/dir/bioserv2/prodgbl07/human/assembly/bios'erv2/firstrun 

/dir/bioserv2/prodgbl07/human/assembly/mir/first_run- 

/dir/bacio2/prodgbl07/human/assembly/bacio/second run 



2 9 . MER.G E : 



We need now to merge all the results of- the assembly runs into one file 
for each. type. For this, we will use the script GuidedMerge. In the 
directory $BASEDIR/human/assembly {or. mouse for mouse), make a directory 
called final. Now we need to make a list of all the con tigs we ' Wo. r iced 
on - i.e. a merge of all of the. list files. For this do: 

sort $BASEDIR/human/cluster/*.list > $BASEDIR/huraan/assembly/ final/ 1 is t . full 

The command for the merge must be run under sh (and not csh or tcsh, 
because, we .need the stdout and stderr in different files), so the command 
is {for the embl file) is something like this: 

sh -c 'GuidedMerge -embl $BASEDIR/hu man./ assembly/ final/list, full: list /dir/bacio 

Because of limitations of interactive sh (only 25 6 characters, on a line), 
and to avoid mistakes, there is a script that runs the command, based 
upon the Task. list file. The commands that should be run are: 

run _m erge_command.pl task_list= SBASEDIR/human/assembly/Task". list contig_list= $ 

run_merge__command.pl task_list= $BASEDIR/human/assembly/Task.list contig_list= $ 

run _raerge_ command. pitas k_list= SBASEDIR/human/assembly/Task . list contig_list= $ 

run_merge_command.pl task_list= $BASEDIR/human/assembly/Task.list contig_list= $ 

run_merge_command.pl task_list= SBASEDIR/human/assembly/Task . list c o n t i g _ 1 is t = ' $ 

You should check the, .errs files of «he merges - they should all be with 
the samenumber of lines, and' the only kind of errors that should appear 
is lines like: 



Couldn't find W703 52 . AA7 3 2 1 8 7 (which means this contig wasn't found 

in any of the output files) 



Finished embl .bacio (which means this file was read until it's end). 
30. HIERARCHY: 

The cview file needs to be processed to make the input of the 

cluster viewer. To do this, you need to run the hierarchy program: 

heirarchy dir= SBASEDIR/ human / assmbl y/ final in = cview. final out = cluster view.fi 
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31. INSTALL DATABASE: 

The following files should be moved to the database directories: 
to the. clusters- directory: 

1. rich^fasta. flipped (should be called there rich_fasta). 

2. cluster_view 

3. clusters. EMBL. 
4 . graph 

If any of the files (especially the last 'two in the human case) are 
larger than 2GB,- they have to be split up. into smaller files. This is 
•done like -'this: 

divide_file.pl in= $BASEDIR/human/assmbly/f inal/embl . final outprefix= $ D A T A B 
divide_file.pl in= $BASEDIR/human/assmbly/final/graph. final outpref ix- $ D A T A 



29. MAKE_INDEX: some of the files that were produced are a part of the output 
of the project. These are the files: 

1. r i c h _ f a s t a . f 1 i p p e d 

2. clusters. EMBL 

3. transcripts 

4. cluster_view 

These files need to be indexed for quick access by GenCore utilities 
(important, for Gen Web viewer). For this you only need a running 
version of GenCore , and then run " make_index ". Each such run 
gives you another file (called "in_file.ind"). 

31. And that's all folks! You are now 'ready 'for the Dry biology' phase! 



Statistical Appendix 

In each stage of the clustering several 
information about this particular versi 
database release). These numbers should 
as a web page in Compugen's internal si 
vary, since the clustering flow varies, 
calculated for version 105, and the way 

r a w d a t a 



est 


size 


Is - 1 


est 






est 


s e q s 


g r e p " 




est 1 w c 


- 1 


est 


bases 


grep - 


v " 


A > " r n a | 


w c 


r n a 


size 


Is - 1 


r n a 






r n a 


s e q s 


grep " 


A. ^ II 


r n a I w c 


- 1 


r n a 


bases 


grep - 


V " 


* > " r n a | 


w c 



numbers can be computed which give 
on of clustering (and this particular 

be generated and formatted for display 
t e •. The pages from each version might 

but here is a summary of the numbers 

they were obtained: 



(take third number minus first number) 
(take third number minus first number) 
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clean data 



est 


umsk 


size 


Is -1 


est 


.umsk 








r n a 


urn s k 


size 


Is -1 


rna 


.umsk 








est 


x o u t 


size 


Is -1 


est 


. x o u t 








•r 'n a 


x o u t 


si z e 


Is - 1 


rna 


.xout 








cln 


e s t s 




g r e p 




est. umsk | wc 


- 1 






c 1 n 


r n a s 




g r ep 


ii /\ ^ it 


. rna. umsk | wc 


-1 






cln 


est 


bases 


g r ep 


- v " 


~ > " est .umsk i 


w c 


( t h ird 


number 


cln 


r n a 


bases 


g r ep 


- v " 


A >" rna. umsk | 
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(third 


number - 



first number) 
first number) 



without abundant genes 



est umsk size 


Is -1 


e-s t . u m s k . o k 












rna umsk size 


Is -1 


rna. umsk. ok. 












est xout size 


Is -1 


est. xout. ok 












rna xout size 


Is -1 


rna. xout. ok 












cln ests 


g r e p 


" A > " est. umsk. ok 


I w c - 1 ■ 










cln r n a s 


g r ep 


" A >" rna. umsk. ok 


| wc -1 










cln est bases 


g r ep 


-v ««">» est. umsk 


. o k 1 wc (third 


number 




first 


number) 


cln rna bases 


gr ep 


-v •» * > ■ rna. umsk 


.ok | wc (third 


n urn b e r 




first 


number) 


first clustering 














total contigs 


w c - 1 


stat 












singletons 


uniq 


- c stat I tail { 


first number of 


last 1 i n 


e ) 






doubletons 


u n i q 


-c stat I tail ( 


first number of 


second 1 


i n e from 


bottom) 


tripletons 


u n i q 


-c stat I tail ( 


first number of 


third 1 i 


n e 


from 


bottom) 


quadrupletons 


u n i q 


-c.stat I tail ( 


first number of 


forth li 


n e 


from 


bottom) 


quintupletons 


uniq 


- c stat 1 taii ( 


first number of 


fifth' 1 i 


n e 


f r o m 


bottom) 


biggest contig 


head 


stat ( 


number in first 


line) 








second contig 


head 


stat ( 


number in. second 


line) 








third contig 


head 


stat ( 


number in third 


line) 








fourth contig 


head 


stat. - . ( 


number in fourth 


line) 








fifth contig 


head 


stat ( 


number in fifth 


line) 








contigs > 50 0 


a w k ' 


$ 1 > 500' stat | 


w c - 1 










contigs > 200 


a w k • 


$ 1 > 2 0 0 • stat | 


w c - 1 











second clustering (for human only) 

Same as first - clustering - use the merged f iles 



clone information 

clusters >1 
biggest cluster 
second biggest 
contig -triplets 
contig-pairs 
standalone ctgs >1 
standalone seqs 
avg contigs/cluster 
without single ctgs 
avg seqs/contig 
without singletons 



wc -1 clone__stats 

head clone — stats (number in first line) 
head clone_stats (number in second line) 

uniq -c clone_stats | tail (first no., second line frombott 
uniq -c clone_stata I tail (first number bottom line) 
awk '$6 = = 1 & & $5 > 1* seq.clu | aw k * (print $2}' | sort | u 
awk '$6 = = 1 & & $ 5 , = = 1 ' seq.clu | wc -1 

t o t a 1 _ c o n t i g s / .[ ( c 1 u s t e r s >. 1 ) + ( s . a . contigs>l ) + (s . a . seqs)] 
[total_contigs - (s.a contigs>l) (s.a. seqs)} / (clusters> 

(cln_ests+cln_rnas without abundant seqs) / total_contigs 
( c 1 n_e s t s + c 1 n_ r n a s without abundant seqs - singletons) / (to 
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Clustering instructions 



http://internal/-raveh/work/projects/cluster/instructions.h 



assembly 

clusteravglen 
without single ctgs 
contigavglen 
without singletons 
.singleton avglen 
%alt splicing 
Avg splice alts. 
Where available 

files sizes 

obtained from "Is -1" on tn « following files: 

est.gb- 

rna.gb 

r i c h _ f as t a 

transcripts 

clusters . EMBL 

keywords 

cluster view 



total 
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